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INTRODUCTION 

Oi  major  consideration  driving  modern  aircraft  gas  turbine  engine 
designs  is  the  desire  to  increase  the  thrust-to-weight  ratio.  This,  in  turn, 
requires  more  highly  loaded  compressor  and  turbine  stages  with  concern  for 
optimum  spacing  between  adjacent  bl.de  passages.  With  this  motivation,  the 
fluid  dynamic  processes  in  typical  turbine  and  compressor  blade  passages  and 
stages  have  become  an  item  of  considerable  importance.  The  aerodynamic 
behavior  of  these  components  can  have  a  major  effect  upon  the  overall  engine 
efficiency.  In  addition,  details  of  the  flow  pattern  can  significantly 
influence  structural  integrity  via  unsteady  forces  produced  by  the  relative 
blade  motion  and,  in  the  case  of  the  turbine,  by  influencing  the  leading  edge 
heat  transfer  rates. 

To  date,  considerable  attention  has  been  focused  upon  single  blade 
passage  flows.  For  example,  Dring,  Joslyn  and  Hardin  (Ref.  1);  Wagner,  Dring 
and  Joslyn  (Ref.  2);  Joslyn  and  Dring  (Ref.  3);  and  Shreeve  and  Neuhoff 
(Ref.  4)  have  performed  detailed  investigations  in  compressor  rotors  and 
stators.  Linear  cascades  have  been  considered  by  Hobbs,  Wagner,  Dannenhoffer 
and  Dring  (Ref.  5),  among  others.  In  regard  to  the  single  passage  flow 
field,  Navier-Stokes  procedures  have  been  applied  to  this  problem  by  several 
investigators  Including  Shamroth,  McDonald  and  Briley  (Ref.  6);  Hah  (Ref.  7); 
Weinberg,  Yang,  McDonald  and  Shamroth  (Ref.  8);  Schafer,  Fruhaufer,  Bauer  and 
Guggolz  (Ref.  9);  and  Dawes  (Ref.  10).  These  efforts  include  both  two-dimen¬ 
sional  and  three-dimensional  analyses  (Refs.  7  and  8).  The  two-dimensional 
analyses  of  Refs.  6,  and  8  through  10  are  based  upon  ADI  solutions  of  the 
Navier-Stokes  equations.  Schafer,  et  al.,  (Ref.  9)  considered  Euler 
solutions  for  flow  through  a  cascade  of  NACA-0012  airfoils,  as  well  as  both 
Navier-Stokes  and  "thin  shear  layer"  solution  for  flow  through  a  cascade  of 
NACA-8410  airfoils.  Although  convergence  comments  were  made  in  regard  to 
Fuller  cascade  solutions  presented  in  this  paper;  for  example,  a  cascade  of 
0012  airfoils  at  zero  degree  stagger  appeared  to  require  approximately  400  to 
600  time  steps  to  converge,  the  number  of  time-steps  to  converge  their  "thin 
shear  layer"  or  Navier-Stokes  calculation  was  not  clear.  Dawes  (Ref.  10) 
also  used  an  ADI  approach  to  solve  the  Navier-Stokes  equations  for  flow  in 
transonic  compressor  cascades;  in  this  case  400-600  time  steps  are  required 
to  converge.  Both  Refs.  9  and  10  considered  compressor  type  cascades  with 
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thin  blades  and  relatively  little  turning.  In  contrast.  Refs.  6  and  8,  which 
were  performed  with  earlier  versions  of  the  code  used  here,  considered  both 
compressor  cascades  and  turbine  cascades.  The  latter,  having  blunter  leading 
edges,  higher  stagger  angles  and  more  turning,  become  a  more  difficult  and 
challenging  simulation  problem.  Calculations  performed  recently  with  this 
code  required  80-120  time  steps  to  convergence,  with  run  time  per  time  step 
for  the  grids  used  being  approximately  0.4  CPU  seconds.  Solutions  were 
obtained  in  30  to  70  CPU  seconds  of  Cray-I  time  using  grids  containing  3500 
points.  In  Ref.  7,  Hah  took  a  different  approach,  using  a  control  volume 
relaxation  approach  to  solving  the  Navier-Stokes  equations.  However,  the 
equations  being  solved  by  a  relaxation  procedure  are  limited  to  a  steady  mean 
flow.  It  may  be  possible  to  extend  such  a  procedure  to  unsteady  flow, 
however,  such  an  extension  would  require  iteration  at  each  time  step,  which 
would  increase  run  time  significantly. 

Although  considerable  information  can  be  gained  by  considering  the  flow 
through  a  single  row  of  blades,  a  more  complex  and  potentially  more  important 
problem  is  that  of  the  flow  through  a  stator-rotor  stage.  Such  an 
investigation  would  include  the  potentially  important  interaction  between 
adjacent  blade  passages.  Information  obtained  from  a  stage  investigation 
could  contain  important  flow  physics  not  addressed  in  a  single  blade  row 
Investigation.  To  date,  considerably  less  effort  has  focused  upon  the  stage 
probLem.  In  regard  to  experimental  efforts,  Detroit  Diesel  Allison  has 
studied  various  aspects  of  unsteady  flow  in  a  compressor  stage  with  the 
primary  emphasis  being  the  time-variant  aerodynamics  of  the  stator  as  a 
function  of  rotor-stator  spacing,  solidity,  etc.  (Refs.  11-14).  Dring, 
Joslyn,  Hardin  and  Wagner  (Ref.  15)  have  considered  the  effect  of  axial 
spacing  upon  the  turbine  stage  flow  field. 

In  regard  to  analysis,  the  rotor-stator  stage  calculation  is 
considerably  more  demanding  than  a  single  blade  row  passage  calculation.  In 
the  stage,  the  flow  is  intrinsically  unsteady  and,  therefore,  transient 
accuracy  is  required.  Furthermore,  the  computational  grid  required  for  a 
stage  is  considerably  more  complex  than  that  for  a  single  passage  analysis. 

To  date,  two  stage  calculation  procedures  which  are  based  upon  solution  of 
the  Navier-Stokes  equations  have  been  developed.  They  are  the  procedure  of 
Rai  (Ref.  16)  and  the  procedure  reported  herein.  Both  efforts  considered  the 
same  basic  problem,  that  of  a  time  dependent,  two-dimensional  Navier-Stokes 
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simulation  of  flow  through  a  stage.  The  major  differences  between  the  two 
were  the  computational  coordinate  systems  and  the  numerical  techniques  used. 
Prior  to  considering  the  details  of  the  present  analysis,  it  is  useful  to 
compare  the  similarities  and  differences  of  the  present  approach  with  that  of 
Ref.  lb  in  a  general  manner. 

In  regard  to  the  computational  grids  used,  Rei  used  a  combination  of 
patched  and  overlaid  grids,  as  shown  in  Fig.  1,  which  is  taken  from  Ref.  16. 
The  blade  geometries  are  shown  in  Fig.  la.  Figure  lb  shows  'O'  grids  in  the 
immediate  vicinity  of  each  blade.  Each  'O'  grid  is  embedded  within  an  'Hf 
grid  where  the  'H'  grids  are  shown  in  Fig.  1c.  The  general  patching  and 
overlay  is  shown  in  Fig.  Id.  Although  this  approach  can  give  a  highly 
resolved  grid  through  the  flow  field,  the  overlay  between  regions  1  and  2, 
and  regions  3  and  4,  requires  a  technique  to  specify  boundary  conditions  on 
the  inner  boundaries  of  zones  3  and  4,  and  the  outer  boundaries  of  zones  1 
and  2.  This  basically  requires  some  method  of  applying  "boundary  conditions" 
on  a  computational  line  internal  to  the  flow  field.  Rai  has  successfully 
approached  this  problem  by  an  iterative  procedure  in  which  outer  zone  1  and 
and  2  boundary  points  are  interpolated  from  zone  3  and  4  results,  and  inner 
zone  3  and  4  boundary  points  are  interpolated  from  zone  1  and  2  results. 
Although  this  approach  has  been  successful,  it  does  require  additional 
iteration  within  the  calculation  which,  in  turn,  increases  computer  run 
time.  In  addition,  the  overlapped  multiple  patched  grids  may  require  more 
grid  points  than  a  single  non-overlaid  grid  for  each  blade.  Finally,  the  use 
of  these  internal  line  boundary  conditions  may  introduce  a  stability  limit  on 
the  allowable  time  step  which  can  be  taken.  Along  the  patched  boundaries  a 
somewhat  similar  approach  is  used  in  which  an  iterative  solution  between 
zones  is  combined  with  conservation  of  fluxes  across  a  patch  line.  A  major 
advantage  of  this  approach  is  that  it  allows  a  relatively  small  gap  between 
rotor  and  stator  to  be  considered. 

In  the  approach  detailed  in  the  present  report,  a  key  priority  was 
development  of  an  efficient  stator-rotor  stage  calculation  procedure  which 
would  require  a  minimum  run  time.  Minimization  of  run  time  was  considered 
critical  for  several  reasons.  First,  even  for  the  simplest  stage 
calculation,  which  is  two-dimensional  and  assumes  equal  pitch  In  rotor  and 
stator,  four  to  eight  cycles  are  expected  to  be  required  to  obtain  a  periodic 
solution.  A  cycle  is  defined  as  the  time  required  for  the  rotor  to  move  a 
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distance  equal  to  one  pitch.  The  goal  of  the  present  effort  was  to  develop  a 
stage  Navier-Stok.es  simulation  technique  which  would  lead  to  periodic 
solutions  for  the  two-dimensional,  equal  pitch  configurations  in  a  modest 
amount  of  computer  run  time  to  allow  stage  calculations  to  be  made  on  a 
regular  basis.  Secondly,  the  eventual  interest  is  in  either  three-dimension¬ 
al  stage  simulations  or  simulations  in  which  rotor  and  stator  have  unequal 
number  of  blades,  i.e.,  unequal  pitch.  The  three-dimensional  simulation  is 
expected  to  require  sixty  times  the  computer  time  of  the  simple  two-dimen¬ 
sional  calculation  and  the  two-dimensional  unequal  blade  calculation  is 
expected  to  require  twenty  to  forty  times  the  computer  run  time  of  the  simple 
two-dimensional  calculation.  If  these  are  to  be  approached  on  anything  but  a 
research  demonstration-type  basis,  the  basic  two-dimensional  equal  pitch  case 
must  be  very  efficient  in  use  of  computer  resources. 

Both  a  single  grid  approach  and  a  patched  grid  approach  were  considered. 
The  trade-off  appeared  to  be  the  single  grid  advantage  of  rapid  run  time 
versus  the  small  gap  spacing  ability  of  the  patched  grid  approach.  Since  one 
important  application  of  the  procedure  would  be  compressor  technology,  and 
since  a  compressor  stage  gap  of  approximately  25%  axial  chord  appears  to  be 
viable  with  the  single  grid  approach,  the  single  grid  approach  was  chosen. 

As  will  be  discussed,  this  has  led  to  a  very  efficient  procedure  which  would 
require  only  2  CPU  minutes  of  CRAY  time  per  cycle  and  10  minutes  to  obtain  a 
periodic  solution.  This  appears,  at  the  present  time,  to  be  significantly 
faster  than  the  run  times  achieved  by  Rai  using  a  patched  grid  and  a  similar 
number  of  grid  points.  The  single  deforming  grid  used  in  the  present  effort 
is  shown  in  Fig.  2.  As  can  be  seen,  this  approach  basicaLLy  combines  a  rigid 
'C'  grid  for  the  leading  stator  and  a  rigid  'H'  grid  for  the  following  rotor 
with  a  deforming  section  between  them.  The  primary  advantage  of  this 
approach  is  the  absence  of  internal  boundaries  such  as  that  which  occurs  in 
overlaid  or  patched  grids,  and  the  consequent  expected  decrease  in  run  time. 
The  main  disadvantage  of  the  single  deforming  grid  approach  lies  in  the 
limitation  as  to  how  close  the  rotor  and  stator  could  be  placed.  Obviously, 
a  small  gap  would  Increase  the  Interaction  in  both  the  rotor  and  stator. 
However,  it  is  estimated  for  practical  modern  compressor  blade  shapes,  a  gap 
of  25%  axial  chord  could  be  calculated  with  this  approach,  and  for  turbine 
blades  a  gap  of  40%  axial  chord  could  be  considered.  These  are  close  to  the 
lower  limits  of  the  gaps  currently  used  in  practice. 
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The  second  difference  in  the  approach  concerns  the  numerical  method. 

The  method  of  Rai  is  a  multiple  grid  approach  for  the  "thin  shear  layer” 
approximation  to  the  Navier-Stokes  equations,  in  which  the  calculation  is 
performed  on  several  grids  which  then  are  patched  together  or  overlaid.  The 
method  used  for  each  individual  grid  is  a  factored,  iterative  implicit 
technique  which  is  an  extension  of  the  noniterative  techniques  of  Ref.  18. 

In  brief,  for  each  time  step,  the  technique  solves  the  equations  on  each 
subgrid  by  an  approximate  factorization  technique.  However,  an  iteration  at 
each  time  step  is  incorporated  within  the  method.  These  additional  iteration 
steps  are  used  primarily  to  match  the  solutions  at  patched  and  overlaid 
boundaries.  Although,  use  of  the  additional  iteration  increases  run  time,  it 
does  have  the  beneficial  effect  of  decreasing  splitting  and  linearization 
errors . 

The  approach  used  in  the  present  effort  solves  the  full  time-dependent 
Navier-Stokes  equations  via  a  technique  in  which  a  single  grid  deforms  to 
accommodate  the  relative  motion.  A  single  double  sweep  LBI  procedure  based 
upon  the  technique  of  Refs.  19  and  20  is  used  to  advance  the  solution  each 
time  step.  Although  an  iteration  could  be  incorporated  to  reduce  the 
splitting  error  and  improve  transient  accuracy,  previous  experience  (Ref.  21) 
indicates  that  at  modest  subsonic  Mach  numbers,  M  0.35,  accurate  time- 
dependent  results  can  be  obtained  for  airfoil  type  problems  in  this  frequency 
regime.  Furthermore,  for  many  aircraft  gas  turbine  applications,  the 
limiting  low  Mach  number  of  0.35  for  transient  accuracy  with  the  time  step 
used  here  falls  within  the  range  of  practical  applications.  The  techniques 
of  the  present  effort  appear  to  be  considerably  more  efficient  in  terms  of 
computer  run  time  for  a  cycle,  and  this  represents  a  major  advantage. 

However,  vis-a-vis  the  patched  grid  approach,  the  technique  is  limited  in 
terms  of  minimum  gap. 


ANALYSIS 


Governing  Equations 

The  equations  used  are  the  full  ensemble-averaged  time-dependent  Navier- 
Stokes  equations  which  can  be  written  in  vector  form  as 

Continuity 


—  +  V-pU  =  0  (1) 

3t 


Momentum 

— U  +  V*(pUU)  =  -VP  +  V*(*+iT)  (2) 

3t 

Energy 


i£.h  +  V 4 (  pUh)  =  -V-CQ+qT)  +  H  +  <J>  +  pc  (3) 

at  Dt 

->•  * 

where  p  Is  density  ,  U  is  velocity,  P  is  pressure,  n  is  the  molecular  stress 
t  * 

tensor,  ir‘  Is  the  turbulent  stress  tensor,  h  is  enthalpy,  Q  is  the  mean 

heat  fLux  vector,  0T  is  the  turbulent  heat  flux  vector,  $  is  the  mean  flow 

dissipation  rate  and  c  is  the  turbulence  energy  dissipation  rate.  If  the 

flow  is  assumed  as  having  a  constant  total  temperature,  the  energy  equation 

is  replaced  by 


Tt  =  T  +  -3 —  =  constant  (4) 

2C 

P 

where  Tt  is  the  stagnation  temperature,  q  is  the  magnitude  of  the  velocity, 
and  Cp  is  the  specific  heat  at  constant  pressure.  In  the  cases  considered 
in  this  work,  the  total  temperature  has  been  assumed  constant.  This 
assumption  was  made  solely  to  reduce  computer  run  time,  and  the  code  has  the 
capability  to  solve  an  energy  equation,  as  has  been  demonstrated  in  Ref.  8. 


It  should  be  noted  that  for  time-dependent  flows  such  as  the  rotor-stator 
flow  field,  the  constant  total  temperature  assumption  Is  an  approximation. 

As  discussed  In  Ref.  21,  the  analysis  of  Ref.  22  shows  the  effect  of  unsteady 
flow  Is  proportional  to  the  Mach  number  squared  and  the  frequency.  In  terms 
of  demonstration  results,  such  as  those  presented  here,  at  the  Mach  number 
and  frequencies  considered,  the  qualitative  nature  of  the  flow  field  should 
be  correct.  The  effect  of  this  total  temperature  constant  assumption  on 
quantitative  results  will  require  further  Investigation.  In  any  case,  future 
efforts  can  easily  Include  an  energy  equation. 

A  number  of  terms  appearing  In  equations  (1-4)  require  definition.  The 
stress  tensor  appearing  in  equation  (4)  is  defined  as 

*  =  2uD  -  u  -  KB)V*UI  (5) 

*  * 
where  KB  is  the  bulk  viscosity  coefficient,  I  is  the  identity  tensor  and  D 

is  the  deformation  tensor,  defined  by: 

D  =  I  ((VU)  +  (VU)T)  (6) 

2 

In  addition,  the  turbulent  stress  tensor  has  been  modeled  using  an 
isotropic  eddy  viscosity  such  that: 


*T 

tr 1  =  - 


pu  u 


$  2  >  t 

—  2ujD  -  —  ( pjV*r) i 


(7) 


where  uj,  the  turbulent  viscocity,  is  determined  by  a  suitable  turbulence 
model.  Turbulence  modeling  is  described  in  the  next  section. 

Equation  (3)  contains  a  mean  heat  flux  vector  defined  as  follows: 


Q  =  -<VT  (8) 


,  AV, 


'  A'V,  .\‘.Y .■  v/v'r 

.  .  iVA  ■  -  «*..  ■■  „•  "."v 
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and  a  turbulent  heat  flux  vector  defined  as: 


qT  =  -<Tvt  (9) 

where  <  and  are  the  mean  and  turbulent  thermal  conductivities. 

Also  appearing  in  equation  (8)  is  the  mean  flow  dissipation  term  4>. 

<t>  =  2uD:D  -(ip-  KB)(V*U)2  (10) 

The  equation  of  state  for  a  perfect  gas 

P  =  pRT  (11) 

where  R  is  the  gas  constant,  the  caloric  equation  of  state 

e  =  CyT  (12) 

and  definition  of  static  enthalpy 

h  =  CpT  (13) 

supplement  the  equations  of  motion. 

In  regard  to  molecular  viscosity,  thermal  conductivity  and  bulk 
viscosity,  the  existing  code  sets  bulk  viscosity  to  zero  and  allows  u  and  < 
to  either  be  constant  or  have  a  Sutherland  type  law  variation.  The  results 
presented  here  were  obtained  with  p  constant  and  the  total  temperature 
assumption  replacing  the  energy  equation. 

Dependent  Variables  and  Coordinate  Transformation 

The  governing  equations,  equations  (1-2),  are  written  in  general  vector 
form  and,  prior  to  their  application  to  specific  problems,  It  is  necessary  to 
decide  upon  both  a  set  of  dependent  variables  and  a  proper  coordinate  trans¬ 
formation.  Based  upon  previous  Investigations  (e.g..  Ref.  b),  the  specific 


scalar  momentum  equations  to  be  solved  are  the  x  and  y  Cartesian  momentum 
equations.  The  dependent  variables  chosen  are  the  physical  Cartesian 
velocities  u  and  v  and  the  density  p. 

The  present  application  focuses  upon  two-dimensional  time-dependent 
flows  with  the  assumption  of  constant  total  temperature.  Therefore,  the 
remainder  of  the  equation  development  will  be  within  this  framework.  The 
equations  are  transformed  from  the  Cartesian  variables  x,  y  to  a  set  of 
general  coordinates  S>  n  where 


S  =  S(x>  y.  t) 

n  =  n(x,  y,  t) 


(14) 


T  =  t 


The  equations  themselves  can  be  expressed  in  either  the  so-called 
'strong  conservation  form" 
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where 


I)  =  5XUy  “  Cy  nx 

(17) 

p  pu  pv  0  0 

W  =  (  pu  )  ,  F  =  (  pu2+p  )  ,  G  =  (  puv  )  Fj  =  (  Txx  )  Gj  =  (  xxy  ) 
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The  advantages  of  the  two  possible  forms  have  been  discussed  in  detail  in 
Refs.  23-25.  If  the  strong  conservation  form  of  the  equations  is  to  be  used, 
then  care  must  be  taken  to  evaluate  the  metric  data  by  a  method  which  is 
consistent  with  a  control  volume  approach  (Refs.  23  and  24).  Usually  this 
requires  numerical  evaluation  of  the  metric  data  even  if  an  analytical 
functional  relationship  for  the  transformation  is  available.  The  analytical 
representation  of  the  metric  data,  £x,  £y»  etc.,  when  combined  with  the 
strong  conservation  form  of  the  equations,  leads  to  significant  error  for  as 
straightforward  a  calculation  as  low  Reynolds  number  flow  about  a  circular 
cylinder  (Ref.  23).  The  quasi-linear  form  of  the  equations  is  much  less 
sensitive  to  the  form  of  metric  evaluation  and  gives  good  results  for  both 
numerical  and  analytical  evaluations  of  the  metric  data.  Furthermore,  good 
results  have  been  obtained  with  this  approach  by  the  present  authors  for  a 
variety  of  steady  state  and  time-dependent  problems  (e.g.,  Refs.  6,  8  and 
21).  Therefore,  the  set  of  equations  represented  by  Eq.  (16)  was  used  in  the 
present  study. 

It  should  be  noted  that  the  previous  calculations  were  made  with  non- 
deforming  coordinate  grids.  In  the  present  stator-rotor  application,  most  of 
the  grid  is  non-deforming,  however,  a  deforming  region  is  present  between  the 
blades.  Therefore,  the  use  of  the  "quasi-linear”  or  "chain  rule 
conservation"  form  with  the  deforming  grid  between  blades  is  an  item  which 
may  require  further  invest  igat ion . 

The  Turbulence  Model 

Although  the  code  contains  both  a  mixing  length  model  and  a  two-equation 
turbulence  model,  and  results  have  been  obtained  with  both  (e.g.,  Ref.  8), 
the  present  effort  utilized  a  mixing  length  type  model.  The  mixing  length 
model  assumes  the  existence  of  a  mixing  length,  and  then  relates  an  eddy 
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viscosity  Mf,  to  the  mixing  length  by 


m 


S' 


!»♦! 


iMt 


vT  -  pe2  (  (^ii  +  iiii)  ) 

dXj  dXi  dXj 


For  flow  regions  upstream  of  the  leading  edge  where  the  flow  is  attached,  the 
mixing  length  is  determined  by  the  usual  boundary  layer  formulation 


l  =  xyD 


^  ^max 


where  <  is  the  von-Karman  constant,  D  is  a  sublayer  damping  factor  and  fc^x 
is  taken  as  0.09  6  where  6  is  the  boundary  layer  thickness.  The  damping 


factor,  D,  which  has  for  the  most  part  been  utilized,  is  the  van  Driest 


damping  factor 


D  =  (1  -  e-y  /27) 


where  y+  is  the  dimensionless  coordinate  normal  to  the  wall,  yuT/v> 

When  the  mixing  length  formulation  is  used  in  a  boundary  layer 
environment,  6  is  usually  taken  as  the  location  where  u/ue  =  0.99. 
However,  this  definition  assumes  the  existence  of  an  outer  portion  of  the 
flow  where  ue  is  independent  of  distance  from  the  wall,  and  assumes  that 
the  location  where  ue  becomes  independent  of  distance  from  the  wall  marks 
the  end  of  the  viscous  region.  In  a  cascade  Navier-Stokes  calculation  no 


such  clear  flow  division  occurs.  Therefore,  the  boundary  layer  thickness,  6, 


is  set  by  first  determining  umax,  the  maximum  velocity  at  each  given 
streamwise  station,  and  then  setting  6  by 


6  =  2.0y(u/Umax  ,  <l} 


i.e.,  <5  is  taken  as  twice  the  distance  from  the  wall  to  the  location  where 


u/ umax  =  ki .  In  the  present  calculation  kj  was  set  to  0.90. 


2 


In  the  wake  region,  the  free  turbulent  shear  flow  model  of  Rudy  and 
Bushnell  (Ref.  26)  is  used.  In  brief,  this  model  evaluated  two  length 
scales,  6X  and  6XX 

5I  =  1*2  x  (y>95  -  yCL)  +  (x  -  xTE)  tan  8°  (22) 

5II  =  2  x  1.2  x  (y. 95  -  ycL) 

where  y.95  is  the  location  where  the  velocity  is  95  percent  of  the  free- 
stream  velocity,  yQE  is  the  centerline  location  and  x^e  is  the  trailing 
edge  location.  This  model  was  developed  for  symmetric  wakes  and  was  modified 
for  the  present  application  by  setting 

(y.95  -  YCL>  =  0.5  x  (y+. 95  -  y'.gs)  (23) 

where  y+#g5  and  y'.gs  represent  the  upper  and  lower  locations  where 

the  wake  veLocitv  reaches  95%  of  the  freestream  velocity.  Following  Ref.  26, 

the  mixing  length,  l,  was  taken  as 

i  -  0.07  x  MIN  (<5t,  %)  (24) 

In  some  preliminary  unpublished  work  at  SRA,  this  model  has  given  good 
agreement  with  data  for  the  near  wake  region  of  a  single  row  cascade.  In  the 
stage  problem  examined  here,  the  wake  of  the  leading  blade  impinges  upon  the 
trailing  blade,  and  it  is  not  clear  how  well  the  details  of  this  can  be 
simulated  with  a  mixing  length  type  model.  However,  in  this  initial  effort, 
a  mixing  Length  model  should  give  reasonable  qualitative  results,  with  a  more 
sophisticated  modeL  being  considered  for  future  efforts. 

Boundary  Conditions 

The  authors'  experience  in  solving  Navier-Stokes  equations  has  indicated 
the  important  role  boundary  conditions  play  in  determining  accurate  solutions 
and  rapid  numerical  convergence.  The  present  approach  follows  the  tack  of 
BrlLey  and  McDonald  (Ref.  27).  In  brief,  this  approach  sets  total  pressure 
and  flow  angle  at  the  upstream  inflow  boundary  and  static  pressure  at  the 
downstream  outflow  boundary.  These  represent  the  physical  boundary 


conditions  for  the  governing  set  of  differential  equations.  In  addition,  the 
density  derivative  is  set  to  zero  at  the  upstream  boundary,  and  second 
derivatives  of  both  velocity  components  are  set  to  zero  at  the  downstream 
boundary.  On  the  blade  surface,  no  slip  conditions  are  applied.  These  are 
applied  with  sublayer  resolution;  wall  functions  are  not  used,  as  the  first 
point  off  the  wall  is  generally  within  the  sublayer.  In  addition,  the  wall 
normal  pressure  gradient  is  set  to  zero.  Finally,  periodic  conditions  are 
applied  in  an  implicit  manner  on  the  periodic  line. 

One  item  of  interest  in  regard  to  boundary  conditions  concerns  the 
downstream  constant  pressure  boundary.  Since  this  is  a  reflecting  boundary 
condition,  it  was  decided  to  explore  possible  non-reflecting  boundary 
conditions.  After  surveying  the  literature,  a  condition  based  upon  the  work 
of  Rudy  and  Strikwerda  (Ref.  28),  and  Engquist  and  Majda  (Ref.  29)  was 
incorporated  into  the  code  upon  option.  Cases  were  run  for  a  single  blade 
row,  however,  this  did  not  appear  to  give  any  advantage  in  convergence 
properties  over  the  constant  pressure  boundary.  In  addition,  the  final 
results  were  nearly  identical  to  those  obtained  with  the  constant  pressure 
boundary  condition. 

Numerical  Procedure 


The  numerical  procedure  used  to  solve  the  governing  equations  is  a 
consistently  split  linearized  block  implicit  (LBI)  scheme  originally 
developed  by  Briley  and  McDonald  (Ref.  19).  A  conceptually  similar  scheme 
has  been  developed  for  two-dimensional  MHD  problems  by  Lindemuth  and  Killeen 
(Ref.  30).  The  procedure  is  discussed  in  detail  in  Refs.  19  and  20.  The 
method  can  be  briefly  outlined  as  follows:  the  governing  equations  are 
replaced  by  an  implicit  finite  difference  approximation,  optionally  a 
backward  difference  or  Crank-Nicolson  scheme.  Terms  involving  nonlinearities 
at  the  implicit  time  level  are  linearized  by  Taylor  expansion  in  time  about 
the  solution  at  the  known  time  level,  and  spatial  difference  approximations 
are  introduced.  The  result  is  a  system  of  multidimensional  coupled  (but 
linear)  difference  equations  for  the  dependent  variables  at  the  unknown  or 
implicit  time  level.  To  solve  these  difference  equations,  the  Douglas-Gunn 
(Ref.  31)  procedure  for  generating  alternating  direction  implicit  (ADI) 
schemes  as  perturbations  of  fundamental  implicit  difference  schemes  is 
Introduced,  in  its  natural  extension,  to  systems  of  coupled  linear  difference 
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equations  having  narrow-block  banded  matrix  structures  which  can  be  solved 
efficiently  by  standard  block-elimination  methods. 

The  method  centers  around  the  use  of  a  formal  linearization  technique 
adapted  for  the  Integration  of  initial-value  problems.  The  linearization 
technique,  which  requires  an  implicit  solution  procedure,  permits  the 
solution  of  coupled  nonlinear  equations  in  one  space  dimension  (to  the 
requisite  degree  of  accuracy)  by  a  one-step  noniterative  scheme.  Since  no 
iteration  is  required  to  compute  the  solution  for  a  single  time  step,  and 
since  only  moderate  effort  is  required  for  solution  of  the  implicit 
difference  equations,  the  method  is  computationally  efficient;  this 
efficiency  is  retained  for  multidimensional  problems  by  using  what  might  be 
termed  block  ADI  techniques.  The  method  is  also  economical  In  terms  of 
computer  storage,  and  in  its  present  form  requires  only  two  time  levels  of 
storage  for  each  dependent  variable.  Furthermore,  the  block  ADI  technique 
reduces  multidimensional  problems  to  sequences  of  calculations  which  are 
one-dimensional,  in  the  sense  that  easily  solved  narrow  block-banded  matrices 
associated  with  one-dimensional  rows  of  grid  points  are  produced. 


Artificial  Dissipation 

The  final  item  to  be  considered  concerns  the  use  of  articifial  dissipa¬ 
tion.  Since  the  calculations  were  at  high  Reynolds  numbers  typical  of  normal 
turbomachinery  applications,  it  was  necessary  to  add  "articifial  dissipation" 
terms  to  suppress  central  difference  spatial  oscillations.  Such  "artificial 
dissipation"  could  be  added  via  the  spatial  differencing  formulation  (e.g., 
one-sided  difference  approximations  for  first  derivatives)  or  by  explicitly 
adding  an  additional  dissipative  type  term.  The  present  authors  favor  the 
latter  approach,  since  when  an  additional  term  is  explicitly  added,  the 
physical  approximation  being  made  is  usually  clearer  than  when  dissipative 
mechanisms  are  contained  within  numerical  truncation  errors,  and  further, 
explicit  addition  of  an  artificial  dissipation  term  allows  greater  controL 
over  the  amount  of  nonphysical  dissipation  being  added.  Obviously,  the  most 
desirable  technique  would  add  only  enough  dissipative  mechanism  to  suppress 
oscillations  without  deteriorating  solution  accuracy.  Various  methods  of 
adding  artificial  dissipation  were  investigated  in  Ref.  32,  and  these  were 
evaluated  in  the  context  of  a  model  one-dimensional  problem  containing  a 
shock  with  a  known  analytic  solution  (one-dimensional  flow  with  heat 
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transfer).  The  methods  considered  included  second-order  dissipation , 
fourth-order  dissipation  and  pressure  dissipation  techniques. 

As  a  result  of  this  investigation,  it  was  concluded  that  a  second-order 
anisotropic  artificial  dissipation  formulation  suppressed  spatial  oscilla¬ 
tions  without  impacting  adversely  on  accuracy,  and  could  be  used  to  capture 
shocks  successfully.  In  this  formulation  the  terms 

|y  (Pn_1dx  3*/ 3x) 

are  added  to  the  governing  equations  where  (J>  =  u,  v  and  p  for  the  x-momentum, 
y-momentum,  and  continuity  equations,  respectively.  The  exponent,  n,  is  zero 
for  the  continuity  equation  and  unity  for  the  momenta  equations.  The  dissi¬ 
pation  coefficient,  dx,  is  determined  as  follows.  The  general  equation  has 
an  x-directlon  convective  term  of  the  form  a 3 d>/  3x  and  an  x-direction 
diffusion  term  of  the  form  3(b34>/ 3x)/ 3x.  The  diffusive  term  is  expanded 

3(b3<f>/3x)/3x=b32<j>/3x2+3b/3x3(j>/  3x  (25) 

and  then  a  local  cell  Reynolds  number  Re^x  is  defined  for  the  x-direction 
by 

ReAx  =  j a- 3b/ 3x| Ax/b  (26) 

where  b  is  the  total  effective  viscosity,  including  both  laminar  and 
turbulent  contributions,  and  Ax  is  the  grid  spacing.  The  dissipation 
coefficient  dx  is  nonnegative  and  is  chosen  as  the  larger  of  zero  and  the 
local  quantity  Ue( axRe^x-l ) .  The  dissipation  parameter  ax  is  a 
specified  constant  and  represents  the  Inverse  of  the  cell  Reynolds  number 
below  which  no  artificial  dissipation  is  added.  The  dissipation  coefficient, 
dy,  is  evaluated  in  an  analogous  manner  and  is  based  on  the  local  cell 
Reynolds  number  Re^y,  grid  spacing  Ay  for  the  y-direction  and  the  specified 
parameter  Oy, 

The  question  naturally  arises  as  to  the  value  to  be  used  for  ox  and 
ov.  Based  upon  a  wide  variety  of  experience  (e.g..  Refs.  6,  8,  21,  32),  it 
has  been  concluded  that  values  of  ax  and  Oy  of  0.05  are  sufficient  to 
suppress  spurious  oscillations,  but  do  not  lead  to  non-physical  damping  of 
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the  mean  flow.  In  some  cases  where  gradients  are  large  and  grids  relatively 
coarse  In  terms  of  the  change  of  dependent  variable  from  point  to  point;  for 
example,  the  leading  edge  stagnation  region;  it  may  be  necessary  to  increase 
the  streamwise  dissipation  to  0.5  locally.  However,  the  normal  dissipation 
is  held  at  0.05. 


RESULTS 


Assessment  of  Time-Accuracy 


The  Navler-Stokes  cascade  code  which  is  used  in  the  present  stator-rotor 
simulation  has  been  successfully  exercised  for  a  variety  of  steady  state 
cascade  flow  fields  (e.g.,  Refs.  6,  8,  23),  however,  an  assessment  was  still 
required  as  to  the  time  accurate  nature  of  the  procedure.  In  this  regard, 
two  calculations  are  relevant.  The  first,  an  oscillating  airfoil  calcula¬ 
tion,  was  performed  by  SRA  under  contract  to  NASA  Langley  Research  Center  as 
part  of  a  dynamic  stall  analysis  (Ref.  21).  The  second  case  considered  was 
performed  under  the  present  effort.  The  case  was  a  25.5%  thick  Joukowski 
airfoil  entering  a  gust  in  which  the  vertical  velocity  is  0.25  times  the 
freestream  velocity. 

In  regard  to  the  oscillating  airfoil  calculation,  the  case  considered 
was  that  of  an  NACA  0012  airfoil  oscillating  sinusoidally  in  pitch  in  a 
stream  with  a  ReynoLds  number  based  on  chord  of  2.08  x  10 6  and  a  Mach  number 
of  0.30.  The  airfoil  oscillated  between  4°  and  20°  with  a  dimensionless 
frequency  of  0.125.  The  case  corresponds  to  Data  RUN  51.005  of  the  data  of 
St.  Hilaire  and  Carta  (Refs.  33  and  34).  Although  this  calculation 
represents  an  isolated  airfoil  rather  than  a  cascade,  It  was  made  with  a  code 
essentially  identical  to  the  cascade  deck  used  here  with  boundary  conditions 
appropriate  to  the  isolated  airfoil  applied. 

Comparisons  between  calculated  and  measured  surface  pressure  coeffi¬ 
cients  are  shown  in  Figs.  3-10.  Three  comparisons  during  the  upstroke  are 
shown  in  Figs.  3  and  4.  As  can  be  seen,  the  agreement  is  good.  The  data  was 
reconstructed  from  the  Fourier  coefficient  given  by  St.  Hilaire  and  Carta 
(NASA  CR-165927).  The  third  measured  data  point  on  the  pressure  surface 
(x/c  “  .066)  gave  very  erratic  results  and  was  not  plotted  for  most  of  the 
comparisons.  The  excellent  comparisons  shown  in  Figs.  3  and  4  give  evidence 
to  the  time-accurate  calculation  for  the  surface  pressure.  Figure  5  presents 
a  comparison  at  a  =  17.7°,  a  >  0.  This  is  near  the  incidence  where  stall 
would  first  be  inferred  from  the  lift  and  moment  curves.  The  figure  shows 
some  discrepancy  between  predicted  and  measured  values  as  the  data  presents 
some  evidence  of  a  vortex  being  shed  on  the  suction  surface  leading  edge. 

The  discrepancy  increases  in  Fig.  6  where  the  data  clearly  Indicates  stall. 
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The  calculated  plateau  on  the  suction  surface,  x/c  “  .15,  seems  to  Indicate  a 
vortex  being  Initiated.  Furthermore,  the  calculated  maximum  suction  peak  at 
a  =  19.5°,  Fig.  6,  Is  considerably  less  than  that  at  a  »  17.7°,  Fig.  5. 

Based  upon  the  plateau  and  the  drop  In  suction  peak,  the  calculated  distribu¬ 
tion  at  19.5°  appears  to  be  beginning  the  stall  process.  The  data  at  19.5° 

Is  presented  with  the  calculation  at  a  =  19.9°,  it  >  0  in  Fig.  7.  Although 
these  are  at  different  values  of  a,  they  represent  pressure  distributions  at 
approximately  the  same  incremental  time  after  stall  is  Initiated;  the 
distributions  are  remarkably  similar.  Comparisons  over  the  downstroke  are 
given  in  Figs.  8-10.  Obviously,  the  basic  trends  are  In  agreement,  as  a 
strong  qualitative  comparison  is  shown  between  the  calculation  and  the 
measured  data.  The  results  of  this  calculation,  particularly  prior  to  stall. 
Figs.  3-5,  indicate  the  time-accuracy  of  the  procedure,  at  least  as  far  as 
surface  pressure  is  concerned.  Although  the  reduced  frequency  for  the 
oscillating  airfoil  calculation  is  approximately  an  order  of  magnitude  less 
than  that  found  in  a  typical  rotor-stator  interaction,  the  accurate  and 
efficient  calculation  of  transients  is  more  demanding  from  a  numerical  point 
of  view  for  low  reduced  frequencies  than  for  high  reduced  frequencies. 

The  second  case  considered  was  performed  under  the  present  effort.  The 
case  was  a  25.5%  thick  Joukowski  airfoil  entering  a  gust  in  which  the 
vertical  velocity  is  0.25  times  the  freestream  velocity.  This  case  was  also 
calculated  by  Giesing  via  an  Inviscid  analysis  (Ref.  35)  and  although  viscous 
effects  are  present  in  the  present  analysis,  it  is  expected  that  the  calcula¬ 
ted  pressure  gradients  should  be  in  reasonable  agreement  if  little  or  no 
separation  occurs.  Calculated  surface  pressure  distributions  from  both  tech¬ 
niques  are  presented  in  Figs.  11-14,  where  T  indicates  the  location  of  the 
gust  relative  to  the  airfoil  leading  edge.  As  can  be  seen,  the  two  calcula¬ 
tions  are  In  good  agreement.  The  major  discrepancies  are  the  appearance  of 
surface  pressure  discontinuities  in  the  inviscid  solution  which  are  not  found 
In  the  viscous  solution,  such  as  at  x/c  =  .25  in  Fig.  12.  However,  this  Is 
the  expected  result  of  viscous  effects.  This  comparison  along  with  the 
dynamic  stall  calculations  of  Ref.  21  serve  to  confirm  the  time-accuracy 
capability  of  the  cascade  code  in  prediction  of  surface  pressure  for  those 
frequencies,  geometries,  Reynolds  numbers  and  Mach  numbers. 


H-Grid  Assessment 

The  stator  system  to  be  used  consists  of  a  ' C'  grid  coordinate  system 
for  the  leading  stator  blade  followed  by  an  'H'  grid  system  for  the  rotor; 
these  two  grids  are  joined  by  a  deforming  section,  as  shown  in  Fig.  2. 
Although  considerable  work  had  been  performed  in  a  ’C’  grid  context,  the 
existing  coordinate  system  had  not  been  exercised  for  an  ' H’  grid.  There¬ 
fore,  prior  to  performing  a  stage  calculation,  it  was  decided  to  perform  a 
single  blade  row  calculation  with  an  *H'  grid  coordinate  system.  The  case 
originally  chosen  was  flow  about  a  25.5%  thick  Joukowski  airfoil  at  zero 
incidence.  This  case  was  chosen  since  the  relatively  thick  leading  edge 
region  presented  a  relatively  difficult  test  case  for  the  'H'  grid.  In 
addition,  the  symmetric  cascade  should  give  a  symmetric  flow  field  which 
would  aid  in  verifying  the  code.  The  'H'  grid  used  is  shown  in  Fig.  15.  A 
comparison  between  pressure  distributions  as  calculated  by  the  ' H'  and  ' C' 
grids  is  shown  in  Fig.  16.  These  were  performed  for  turbulent  flow  and 
showed  very  good  correspondence.  The  'C'  grid  capability  has  been  assessed 
in  a  variety  of  comparisons  (e.g..  Refs.  6,  8,  21,  23)  and  therefore,  this 
serves  to  confirm  the  'H'  grid  capability.  In  particular,  it  should  be  noted 
that  the  stagnation  point  pressure  coefficient  with  the  'H'  grid  was 
approximately  0.96.  This  Is  at  the  coordinate  singularity,  and  represents  a 
demanding  computation  for  the  'H'  grid.  The  computed  Cp  of  0.96  should  be 
regarded  as  quite  good.  In  addition,  the  calculation  proved  to  be  symmetric, 
as  indicated  by  the  contour  plot  shown  in  Fig.  17. 

Stator-Rotor  Calculation 

The  final  case  of  interest  is  the  stator-rotor  stage  calculation,  which 
is  the  focus  of  the  present  effort.  The  case  chosen  for  calculation  was  a 
stage  formed  from  two  sets  of  compressor  blade  profiles  of  Hobbs,  Wagner, 
Dannenhoffer  and  Dring  (Ref.  5).  This  blade  shape  has  been  previously  con¬ 
sidered  as  a  cascade  at  SRA,  and  numerical  simulations  of  the  flow  in  this 
cascade  had  shown  good  correspondence  when  compared  against  the  experimental 
data  (Ref.  6).  The  calculation  was  run  for  a  stag'  in  which  the  Reynolds 
number  at  inflow  was  6.R  x  10^;  the  Inflow  Mach  number  was  0.47.  The  grid 
used  consisted  of  a  constructive  ’C'-grld  for  the  upstream  stator,  a 
constructive  ’H’  grid  for  the  downstream  rotor  and  a  taut  spline  blending 
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region  between  the  two.  As  previously  discussed,  the  calculation  set 

upstream  total  pressure  and  flow  angle  and  downstream  static  pressure.  The 

upstream  flow  angle  was  set  at  a  value  of  38°  and  the  upstream  total  to 

downstream  static  pressure  set  at  a  ratio  of  1.0807.  The  rotor  wheel  speed 

was  set  at  52  per  cent  of  the  velocity  upstream  of  the  stator  leading  to  the 
incidence  on  the  rotor  in  a  frame  tied  to  the  rotor  being  approximately  38°. 

A  figure  showing  these  velocities  is  presented  in  Fig.  18. 

The  calculation  was  run  on  a  very  highly  resolved  grid.  The  first  point 
off  the  blade  was  approximately  3xl0-5  chords  from  the  blade.  High 
resolution  was  obtained  in  both  the  stator  and  rotor  leading  edge  regions. 

The  blade  spacing  between  rotor  and  stator  was  0.5  axial  chords;  the  spacing 
between  adjacent  rotor  blades  (or  stator  blades)  was  0.72  axial  chords.  A 
computer  plot  of  the  grid  is  shown  at  three  different  times  in  Fig.  19.  As 
can  be  seen,  the  grid  is  a  rigid  '  C'  grid  in  the  vicinity  of  the  stator; 
i.e.  upstream  of  point  A,  and  a  rigid  'H'  grid  in  the  vicinity  of  the  rotor; 
i.e.  downstream  of  point  B.  Between  A  and  B  the  grid  deforms.  A  cycle 
begins  with  the  grid  in  the  position  of  the  lower  figure,  T  =  8.0.  The 
calculation  proceeds  with  the  grid  moving  until  the  position  at  T  =  9.0  is 
reached.  This  represents  a  cycle;  the  grid  is  then  moved  back  to  the 
position  at  T  =  8.0,  the  results  in  the  deformed  region  are  interpolated,  and 
the  calculation  then  continues. 

The  calculation  required  approximately  five  cycles  to  reach  a  condition 
of  cyclic  periodicity  in  the  solution.  Each  cycle  required  100  time  steps 
and  approximately  4  CPU  minutes  of  CRAY  run  time.  The  code  being  used  is  a 
partially  vectorized  code;  a  more  recent  version  of  this  same  code  is  fully 
vectorized  and  runs  twice  as  fast  as  that  used  here.  Therefore,  it  is 
reasonable  to  expect  run  times  of  2  CPU  CRAY  minutes  per  cycle  for  7500  grid 
points  to  he  attainable.  Plots  of  Mach  number  contours  for  this 
configuration  at  three  different  time  levels  are  shown  in  Figs.  20  and  21. 
During  this  calculation,  minor  changes  were  made  in  the  turbulence  model,  as 
well  as  In  the  leading  edge  region  artificial  dissipation.  Therefore,  this 
does  not  represent  nine  cycles  of  calculation.  As  previously  stated,  a 
reasonably  periodic  solution  is  obtained  within  five  cycles.  For  the  current 
test  case,  the  major  interactive  effect  is  on  the  trailing  blade.  This  can 
be  seen  clearly  in  the  Mach  number  contours  where  the  flow  field  changes 
significantly  between  T  =  N  and  T  =  N+0.5.  In  the  latter  case,  the  stator 
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wake  Impinges  directly  upon  the  rotor  leading  edge.  In  the  former  case,  it 
lies  within  the  passages.  No  major  effect  is  noted  on  the  flow  in  the 
upstream  blade  passage  for  this  blade  streamwise  spacing.  The  change  in 
static  pressure  contours  through  the  cycle  appears  less  dramatic  than  the 
change  in  Mach  number  coutours.  However,  significant  changes  do  occur 
particularly  in  the  leading  edge  regions  when  viewed  in  terms  of  pressure 
coefficient  perturbation. 

A  plot  of  the  'baseline'  rotor  surface  pressure  coefficient  is  shown  in 
Fig.  22.  This  represents  the  surface  pressure  at  a  time  when  the  rotor 
leading  edge  location  is  midway  between  the  stator  trailing  edges.  It  should 
be  noted  that  this  should  be  approximately  the  minimum  interference 
configuration.  It  does  not  represent  a  result  time-averaged  over  the  cycle. 
The  rotor  leading  edge  stagnation  point  pressure  coefficient  is  considerably 
greater  than  unity  (“1.33)  due  to  the  work  being  performed  on  the  fluid  by 
the  moving  blade.  The  pressure  coefficient  was  smooth  throughout  the  flow 
except  in  the  immediate  vicinity  of  the  leading  edge  (x/c<.001)  where  some 
oscillations  occurred.  In  any  case,  it  is  confined  to  a  very  small  region 
of  the  flow. 

Pressure  coefficient  perturbations  relative  to  the  base  case  are 
presented  in  Figs.  22-28.  In  these  figures,  the  relative  positions  of  the 
rotor  leading  edge  and  stator  trailing  edge  art  noted  on  the  left  side  of  the 
figures.  It  should  be  noted  that  surface  pressure  perturbations  are  expected 
to  occur  as  the  rotor  passes  through  the  nonuniform  exit  flow  field  generated 
by  the  stator.  As  the  rotor  passes  through  the  stator  wake  it  experiences 
(i)  a  decrease  in  flow  stagnation  pressure  due  to  losses  in  the  wake  and 
(ti)  an  Increase  in  flow  incidence  due  to  a  decrease  in  axial  velocity.  Tt 
is  estimated  that  the  flow  Incidence  changes  from  approximately  32°  from 
horizontal  when  the  rotor  is  between  stator  blades,  to  approximately  38°  when 
the  rotor  is  in  the  stator  wake.  In  addition,  when  the  rotor  is  in  the 
stator  wake,  it  experiences  an  oncoming  shear  flow  which  will  affect  the 
blade  loading. 

The  results  of  Figs.  23-28  show  the  surface  pressure  perturbation  as  the 
rotor  passes  through  the  cycle.  The  base  case  shown  in  Fig.  22  gives  the 
pressure  distribution  when  the  rotor  is  between  the  stator  blades.  At  this 
location,  the  incoming  flow  angle  relative  to  the  rotor  is  at  its  minimum 
angle  to  the  horizontal.  As  time  progresses,  the  rotor  moves  into  the  wake 


and  the  relative  flow  angle  Increases,  giving  an  Increased  pressure  on  the 
leading  edge  region  pressure  surface  and  a  decreased  pressure  on  the  leading 
edge  region  suction  surface.  This  becomes  most  pronounced  at  T  =  8.59  to 
8.67,  as  shown  In  Figs.  26  and  27.  During  this  portion  of  the  cycle 
considerable  additional  leading  edge  loading  occurs  with  ACp  being  approxi¬ 
mately  0.'.  At  all  times  the  major  loading  perturbations  are  confined  to  the 
leading  edge  region,  X/C< 0.2.  As  the  flow  proceeds  into  the  passage, 

X/OO.2,  the  flow  basically  behaves  as  a  channel  flow  and  the  surface 
pressure  distribution  Is  relatively  Insensitive  to  the  incoming  flow  angle. 
The  major  effect  being  In  the  leading  edge  region  has  serious  Implications  In 
terms  of  time-dependent  pitching  moment  and  possible  structural  consequences. 

Although  the  major  effect  is  in  the  rotor  leading  edge,  there  is  a 
smaller  effect  over  the  entire  blade.  As  shown  in  Figs.  23-28,  both  the 
pressure  and  suction  surfaces  show  a  pressure  perturbation  over  the  entire 
bLade  chord.  Since  this  perturbation  becomes  very  small  at  T  »  9.0,  Fig.  28, 
it  does  not  represent  a  significant  aperlodiclty  still  present  In  the 
solution.  Downstream  of  the  leading  edge  the  perturbation  Is  more  pronounced 
on  the  suction  surface,  reaching  a  maximum  perturbation  Cp  of  -0.025. 

The  perturbation  flow  field  in  the  trailing  rotor  blade  was  considerably 
stronger  than  that  observed  on  the  stator.  However,  some  relatively  small 
effects  were  noted  on  the  stator.  For  example,  the  maximum  suction  peak  on 
the  stator  ranged  from  -0.86  to  -0.81,  giving  a  ACp  of  approximately  0.05, 
with  the  highest  value,  -0.81,  occurring  when  the  trailing  blade  was 
inprox i mate  ly  aligned  with  the  Leading  blade.  The  trailing  edge  pressure 
coefficient  ranged  from  0.34  to  0.37,  with  the  highest  values  again  being 
when  the  stator  trailing  edge  and  rotor  leading  edge  were  approximately 
aligned.  These  results  indicate  a  maximum  interaction  effect  occurs  on  the 
leading  stator  blade  when  the  blades  are  nearly  aligned,  and  minimum 
interaction  occurs  when  the  rotor  is  between  stators.  However,  it  should  be 
noted  that  these  stator  perturbations  are  relatively  small. 

In  regard  to  data  comparison,  no  data  exists  for  this  specific  case. 
However,  the  data  of  Fleeter,  Jay  and  Bennett  (Ref.  36)  can  be  used  for 
guidance.  This  data  was  taken  on  a  rotor/stator  configuration  with  the  rotor 
being  the  leading  blade.  The  gap  was  approximately  0.4  axial  chords.  In 
general,  this  data  shows  the  major  perturbation,  Cp,  to  be  confined  to  the 
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leading  edge  region  of  the  trailing  blade.  When  normalized  by  the  incoming 
dynamic  head,  the  leading  edge  perturbation  was  approximately  0.10  to  0.30 
over  a  range  of  flow  conditions.  This  is  the  same  order  found  in  the  present 
simulations.  As  in  the  present  numerical  flow  simulation,  although  some 
perturbation  was  found  downstream  of  the  leading  edge,  it  was  considerably 
less  than  at  the  leading  edge.  This  again  is  consistent  with  the  results 
obtained  here. 


CONCLUSIONS  AND  RECOMMENDATIONS 


the  present  effort  has  led  to  the  development  of  a  numerical  procedure 
of  the  stage  flow  field  based  upon  solution  of  the  Navler-Stokes  equations. 
The  specific  case  considered  is  the  basic  case  In  which  two-dimensional  flow 
and  equal  rotor  and  stator  pitch  are  assumed.  The  procedure  developed  is 
capable  of  obtaining  periodic  solutions  for  a  grid  of  7500  points  within  ten 
CPU  minutes  of  Cray  run  time.  This  represents  a  very  efficient  technique 
which  will  permit  runs  of  this  type  on  a  regular  basis,  and  which  will  allow 
extension  to  either  three  dimensions  or  to  stages  in  which  rotor  and  stator 
have  unequal  pitch  and  to  multiple  stages. 

The  results  obtained  show  the  qualitatively  expected  features. 
Periodicity  was  obtained  essentially  within  five  cycles.  The  major  effect  of 
the  interaction  to  the  specific  case  considered  appeared  in  the  leading  edge 
region  of  the  downstream  rotor  blade.  Significant  perturbation  pressures 
were  noted,  ACp^O-A,  which  has  major  implications  regarding  unsteady  loading 
and  unsteady  fluid  structure  interaction.  The  concentration  of  the 
interaction  effect  in  the  leading  edge  has  particular  significance  in  regard 
to  unsteady  pitching  moment.  Although  less  significant,  unsteady  effects 
were  noted  over  the  aft  porti  n  of  the  trailing  blade,  as  well  as  over  the 
entire  Leading  blade.  -*■ 

Possible  extensions  of  this  work  fall  into  three  natural  paths.  Under 
the  first  additional  two-dimensional,  equal  pitch  calculations  would  be  made 
and  the  results  investigated  with  particular  emphasis  on  the  flow  physics  and 
possible  structural  Implications  of  unsteady  loading.  The  second  path  would 
consider  three-dimensional  flow  in  a  stage.  This  would  allow  inclusion  of 
the  various  vortices  generated  in  the  three-dimensional  passage.  Since  the 
three-dimensional  calculation  is  expected  to  require  sixty  times  the  computer 
resources  as  the  two-dimensional,  an  efficient  two-dimensional  analysis  is 
Imperative  as  a  base.  The  analysis  developed  under  the  present  effort  is 
very  efficient  and  therefore  would  form  a  suitable  base  In  meeting  these 
goals.  The  third  option  is  extension  to  the  case  of  unequal  rotor  and  stator 
pitch.  Again,  for  this  to  be  practical,  an  efficient  base  deck  is  key. 
Finally,  a  more  long-term  goal  would  be  development  of  a  two-dimensional 
multi-stage  analysis  based  upon  this  approach. 
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Pressure  coefficient  comparison  —  NACA  0012  airfoil. 
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FIGURE  f2.  PRESSURE  COEFFICIENT  FOR  0.255  THICK  JOUKOWSKI 
AIRFOIL  ENTERING  GUST,  T  =  0.25 
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Figure  15.  H-Grid  Joukowski  Cascade. 


SSfl 


i 


Figure  16.  Comparison  of  pressure  distributions, 
Joukowski  Cascade. 


Figure  18.  Inflow  and  blade  velocity. 
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